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We provide a band structure with low-energy properties consistent with recent photoemission and 
quantum oscillations measurements on FeSe, assuming mean-field like s and/or d-wave orbital order¬ 
ing at the structural transition. We show how the resulting model provides a consistent explanation 
of the temperature dependence of the measured Knight shift and the spin-relaxation rate. Fur¬ 
thermore, the superconducting gap structure obtained from spin fluctuation theory exhibits nodes 
on the electron pockets, consistent with the ’V’-shaped density of states obtained by tunneling 
spectroscopy on this material, and the temperature dependence of the London penetration depth. 
Our studies prove that the recent experimental observations of the electronic properties of FeSe are 
consistent with orbital order, but leave open the microscopic origin of the unusual band structure 
of this material. 
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The electronic properties and the nature of the inter¬ 
actions that drive the low-energy physics and the ordered 
phases of iron-based superconductors (FeSC) continue to 
pose an outstanding problem in modern condensed mat¬ 
ter physics. The diversity of the properties among the 
different families of FeSC and their complex multi-orbital 
band structure have hindered the understanding of the 
electronic states in these materials, as well as the mech¬ 
anism of superconductivity. 

A material that stands out is the structurally sim¬ 
plest compound, FeSe, which exhibits a tetragonal to or¬ 
thorhombic structural phase transition at Tg ^ 90 K 
without concomitant spin density wave (SDW) order, 
and becomes superconducting below ~ 9K. Below 
Ts, the material exhibits strong electronic anisotropy and 
the absence of tetragonal symmetry-breaking SDW or¬ 
der, makes FeSe ideal for studying the origin and con¬ 
sequences of nematicity, i.e. the breaking of rotational 
symmetry while preserving translational symmetry. For 
example, an early scanning tunneling microscopy (STM) 
study of FeSe films on SiC substrate found highly elon¬ 
gated vortices and impurity states, and an associated 
nodal superconducting gap,[T] but until recently simi¬ 
lar experiments on crystals were hampered by sample 
quality. Other remarkable properties of FeSe include the 
significant enhancement of the superconducting critical 
temperature both under pressure, [2] and for monolay¬ 
ers of FeSe grown on SrTiOs surfaces. [3H5] 

Recently, the study of bulk FeSe crystals has been revi¬ 
talized by the growth of very clean samples amenable to 
study the details of the low-energy properties by e.g. nu¬ 
clear magnetic resonance (NMR), transport, STM, angu¬ 
lar resolved photoemission spectroscopy (ARPES), and 
quantum oscillation (QO) experiments.!^ Even though 
a consensus on the electronic bands has not yet been 
reached by ARPES. [TUI^ recent studies found that the 
Fermi surface (FS) above Ts consists of two small hole 


cylinders of mainly dxz/dyz character around the V — Z 
line. The hole bands are split by a sizable spin-orbit cou¬ 
pling (SO) of Ape — 20 meV above the structural transi¬ 
tion, and turns into a single elongated hole cylinder below 
TsM ARPES also finds an electron pocket at the M 
point of mainly dxz/dyz character. Importantly, the ex¬ 
pected dxzjdy^ degeneracy at M is lifted by ~ 50meV, 
constituting strong evidence for orbital order in FeSe. 
We emphasize that these results for the electronic struc¬ 
ture are very different from those obtained within DFT 
calculations. Pin] For example, ARPES finds that the 
electronic bands in FeSe are renormalized compared to 
DFT calculations by a factor of ^3 for the dxz/dy^ bands 
and ~9 for the dxy band. pun] QO performed at low T 
in magnetic fields large enough to suppress superconduc¬ 
tivity are consistent with the ARPES data in observing 
small largely 2D pockets, even though the amount of dis¬ 
persion along kz remains unsettled. nnun] 

Recent ^^Se NMR measurements on FeSe have re¬ 
ported a clear splitting of the NMR line shape setting in 



FIG. 1. (Color online) (a) Fermi surface (FS) above Ts with 
SO coupling, and (b) FS at low T with additional orbital 
splitting of 50 meV consisting of a F centered FS (FSrz) cylin¬ 
der and an inner and outer FS centered around the M point 
(FSm.i and FSm,o). The inset shows the T dependence of the 
orbital splitting Aa(T); the two colored (gray) dots represent 
the T chosen for the two displayed FS. 
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at Tg, with an order parameter-like T dependence below 
Ts.[l6] At high r, the spin-lattice relaxation rate is, how¬ 
ever, unaffected by the structural transition, and only ex¬ 
hibits a clear upturn at low T closer to the superconduct¬ 
ing T^. |16L[T7] These recent experiments have been inter¬ 
preted as evidence for orbitally driven nematic behavior 
in FeSe, but despite the apparent weakness of (momen¬ 
tum summed) spin fluctuations near Tg, the spin-nematic 
picture may still apply. One possibility is that, unlike 
other FeSC, fluctuations in FeSe with different wave vec¬ 
tors compete to frustrate long-range magnetic order. [IB] 

Finally, we note that the resulting Fermi energy of 
the bands of FeSe seen by ARPES are remarkably 
small, comparable to the superconducting gap, which 
suggests the possibility that FeSe may be close to a 
BEC/BCS crossover, and thus exhibit unusual thermo¬ 
dynamic properties and magnetic field effects. [TB] Thus, 
for multiple reasons, it is important to perform new the¬ 
oretical studies of this intriguing material and obtain a 
minimal model capturing its main electronic properties. 

Here, we perform a theoretical study of the conse¬ 
quences of orbital order in a band relevant to FeSe. 
Starting from the DFT-generated band for FeSe obtained 
by Eschrig et al. we apply band renormalization of 
Htb = Hq jz (where z = 6 is renormalization factor and 
iJo is the unrenormalized tight-binding Hamiltonian) and 
additional shifts to the hopping integrals (see Supple¬ 
mentary Information (SI) for details) to generate a new 
tight-binding model. We find that a band consistent with 
ARPES and QO is only possible in the presence of an or¬ 
bital splitting setting in at Tg, and a T independent SO 
coupling. Further, we explain the recent Knight shift and 
the spin relaxation rate measurements, and study how 
spin fluctuation-mediated pairing can lead to a nodal gap 
structure in agreement with measured density of states 
(DOS) and penetration depth A of FeSe. 

The bare Hamiltonian used in this study is given by 

H = Htb+Hoo, ( 1 ) 

Htb = V(k)c]),,(k)c,^(k), (2) 

Hqo = As(r) y^(na;zer(k) - Hyzajk)). (3) 

kcT 

Here (/r, J^) are orbital indices, are the hopping 

integrals, and Uf^aik) = cjj^(k)c^o'(k). All details of 
the hopping integrals are provided in the SI for both 
a 5-orbital and 10-orbital model. In the orbitally or¬ 
dered state, Hqo contributes and As(T) is assumed to 
exhibit a mean-field T dependence with a maximum am¬ 
plitude As(T = 0) = 50meV. For simplicity we focus 
in the main part of this paper on a pure s-wave 00, 
but the consequences of an additional d-wave 00 of the 
form Ad{T) Y.k^(cos{k^) - cos{ky)){n^^^{k) + ny^^{k)) 
have also been studied and the results can be found in 
the SI. It has been reported by ARPES that the band 


splitting of the dxz/dyz bands at the M-point in the or¬ 
bitally ordered state does indeed show a mean-field be¬ 
havior and saturates at low T with a band splitting of 
~ 50 meV. [71 [TTl ITB] . Finally, we have included a SO 
term, Hso = Ape J2x,y,z which causes a band 

splitting of 20 meV in the tetragonal high T phase. EHIIS] 
Band Structure. As shown in Fig. and Fig. the 
band structure and the resulting FS of H is in nearly 
quantitative agreement with experiments. Below Tg, the 
hole band at the F point is split by a 50 meV orbital order 
(at T = lOK) and the bottom of the band lies ^ 20meV 
below the chemical potential. Similarly a dispersionless 
dxy-hanA is present at an energy of ~ — 50meV at the 
F point. At the M point, the electron pockets consist 
of quasi-2D cylinders where the outer pocket, having a 
dominant dxy character, encloses an inner dxz/dyz elec¬ 
tron pocket. The inner electron band at the M point 
has an orbital splitting of 50meV and almost grazes the 
Fermi level. These low T band structure values are in 
good agreement with ARPES results. [71 I MTTl ITB] 

Similar agreement with ARPES is achieved for T > Tg 
where the orbital order is absent. There, the hole pock¬ 
ets consist of a quasi-2D outer circular cylinder and an 
inner hole pocket near the Z-point as seen in Fig. [^a). 
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FIG. 2. (Color online) (a) Band structure of the 10 orbital 
model at low T with orbital order and SO coupling, yielding 
the QO frequencies as a function of magnetic held angle 6 
shown in (b), where the error bars indicate the numerical 
uncertainty in the determination of the extremal orbits. The 
orbital character in the band plot is indicated by the colors 
red dxz, green dyz, blue dxy, yellow dj, 2 _j, 2 , purple d^,^ 2 _x 2 . 
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The band also exhibits overall agreement with the orbital 
content observed in polarized ARPES experiments. [TO] 
At high T, the hole pocket at the zone center and the 
inner electron pocket at the M point contain both d^z 
and dyz character, and the outer electron pocket is pre¬ 
dominantly of dxy character. Similar orbital content 
of the Fermi pockets have also been seen in ARPES 
measurements. [m At low T, the orbital content of the 
hole cylinder is dominated by d^z character (dy^ for the 
other twin). For pure s-wave orbital order, both the hole 
and the inner electron cylinders have dominant dxz or¬ 
bital character (see Fig.j^b)) whereas in the presence of 
an additional d-wave orbital order the electron cylinder 
can have the opposite dyz orbital character (see SI). For 
the electron pocket at low T, the inner pocket contains 
both dxz and dyz orbital character whereas the outer 
pocket at low T has orbital content dominated by the 
dxy-orhital. Although the orbital content of the electron 
pockets agree well with experiments. [TOHI^ the outer dxy 
electron pocket is difficult to observe in ARPES due to 
matrix elements effects. 

Quantum Oscillations. The extremal FS areas in FeSe 
at low T as well as their kz dispersion have been stud¬ 
ied by QO measurements. unmang These experiments 
have found four well separated QO frequencies, and ar¬ 
guments have been put forward that the QO frequen¬ 
cies correspond to one electron and one hole quasi-2D 
FS cylinder, m as well the possibility of a single quasi- 
2D hole cylinder and two almost dispersionless electron 
cylinders. [To] Although the former possibility cannot be 
ruled out, in this study we have pursued the latter pos¬ 
sibility, which is supported by the weak dispersion 
observed for the electron cylinders. m 

Starting from the 10 orbital tight-binding Hamiltonian 
and including the effects of SO coupling, we calculate 
the eigenenergies ^i(k) on a grid in the Brillouin zone 
(BZ) and obtain the extremal areas F of the FS for cuts 
on planes perpendicular to the external magnetic field 
using a numerical method. [231 HTj The direction of the 
magnetic field is then parametrized by the angle 9 be¬ 
tween the crystallographic c axis and the field direction. 
For 6 = 0, the electron pockets have extremal areas of 
F ~ 588 T for the dxy pocket and F ^ 102 T for the 
smaller dxz/dyz pocket, as seen from Fig. [^b). The hole 
Fermi cylinder is elongated due to the effect of orbital 
ordering with a maximum area of F ^ 691 T for k^ = tt 
and a minimum area of F ~ 260 T at = 0. Overall the 
experimentally observed binding energies for the dxz/dyz 
and dxy bands, the 3D FS structure of both hole and elec¬ 
tron pockets, the extremal orbit areas as well as their kz 
dispersion are in good agreement with our calculations. 
We have also calculated the Sommerfeld coefficient from 
the effective masses extracted from our quantum oscilla¬ 
tion calculation. Using the recipe given in Ref. [THl we 
find a Sommerfeld coefficient of 4.5 mJ/mol-K^ which is 
in reasonable agreement with experimental value of ^5.3- 


5.7 mJ/mol-K^. [T71 1^ 

Nuclear Magnetic Resonance. Next, we test our elec¬ 
tronic model of FeSe to see if it can also reproduce 
NMR experiments. [m Hr] For computational simplic¬ 
ity, in the following we apply the 5 orbital model that 
shows a good agreement with the 10 orbital model (see 
SI) and ignore the effect of SO coupling which causes 
only small quantitative changes to the observables dis¬ 
cussed in the remainder of this paper. The NMR Knight 
shift is proportional to the homogeneous susceptibility, 
K = AhfXB.PA{fi = 0) -I- ATchem, where we have ap¬ 
proximated the spin susceptibility by its standard RPA 
form, Ahf is the hyperfine form factor, and ATchem is 
a T independent chemical shift which we have ignored 
for the purposes of this study. For the following calcu¬ 
lations, we include local Coulomb interactions via the 
standard Hubbard-Hund Hamiltonian|27j parametrized 
by the Hubbard interaction Uz (z is the band renormal¬ 
ization) and the Hund’s exchange J (see SI), calculate 
the orbitally resolved non-interacting susceptibility, and 
include interactions within RPA.^S] 

In the paramagnetic state, the form factor is a di¬ 
agonal matrix with components (A^^, A)(j, A))j) where 
the coordinates point along the Fe-Fe direction repre¬ 
senting the magnetic field orientation of the NMR ex¬ 
periment. The form factor maintains the symmetry of 
the underlying lattice such that for the high T tetrag¬ 
onal phase A^j = A™ ^ A^j whereas the orbital or¬ 
dered orthorhombic phase has A^^ ^ A//j ^ A^y. This 
anisotropy leads to a split Knight shift frequency below 
Ts in twinned samples of FeSe. [niiiT] 

As shown by Baek et a7,[T5] the Knight shift split¬ 
ting exhibits a T dependence proportional to the mean- 
field orbital order parameter. Therefore, we model the 
form factor by the expression A^f = a ± g{T), where 
g{T) = /3As(T), {a, (3) are fitting parameters, and ± 
refers to the two orthorhombic domains h and I 2 . The 
calculated Knight shift as a function of T is shown in 
Fig- Ha). At high T above Ts, the Knight shift increases 
with T similar to experiments, in contrast to the DFT 
generated non-renormalized bands (see SI). Below Ts, 
for a particular magnetic field direction the Knight shift 
shows a minimum value around T ~ 60 K similar to ex¬ 
perimental results. Below T ~ 60K, we find a slight 
enhancement of the Knight shift signal. Although the 
measured Knight shift saturates and does not show this 
enhancement for both orthorhombic domains, this may 
be simply related to a deviation of the splitting from 
mean-field behavior found experimentally at the lowest 

T.m 

In order to study the evolution of the spin fluctuations. 
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we have also calculated the spin-lattice relaxation rate 
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NMR experiments probing the ^"^Se atoms in FeSe ex¬ 
hibit a q dependent hyperfine form factor in the param¬ 
agnetic state given by Aj^{q) = cos{qx/2) cos(gy/2) 
assuming that the Se ion interacts with its four near¬ 
est Fe neighbors only. Since ^^Se is a spin 1/2 ion, 
quadrupole type coupling to local lattice distortions do 
not contribute to the relaxation rate. As seen from the 
form factor, spin fluctuations at the edges of the BZ will 
be filtered out. The result of the calculation for 1/TiT 
for interaction parameters Uz = 1.8 eV and Jz = QAUz 
is shown in Fig. [^b). As seen, the spin fluctuations are 
enhanced at low T. However, as observed in recent NMR 
experiments [T3 I17j the enhancement does not occur at 
Ts despite the sharp increase of Ag at Tg, but below 
about T ~ 40K. Interestingly, this increase of spin fluc¬ 
tuations at low T is caused by the orbital ordering which 
leads to a low-T incommensurability in the spin suscep¬ 
tibility that pushes spectral weight away from BZ edges, 
and therefore does not allow the structure factor to effec¬ 
tively filter out those fluctuations (see SI). Note that al¬ 
though the low T spin susceptibility avoids the magnetic 
state by remaining below the Stoner limit, the enhanced 
fluctuations at low T have important consequences for 



T(K) 


FIG. 3. (Color online) (a) NMR Knight shift versus T. The 
hyperfine form factor has been taken as = 0.6[0.57 ± 

0.035(/(r)] for 77 T c and A/j. = 0.6 x 0.5 for 77||c. (b) Spin- 
lattice relaxation rate versus T with A/Y^^ = 0.57±0.035(;(r) 
for 77 T c and A^f = 0.5 for 77||c. Red curve H c (domain 
h), Green curve 77 T c (domain I 2 ), Black curve H A c 
(domain average). Blue curve 77||c. 




FIG. 4. (Color online) Superconducting order parameter as 
calculated from spin-fluctuation pairing using the interactions 
Uz = 1.8 eV and Jz = Q.lUz shows nodal regions on one 
of the electron pockets (a). The corresponding DOS clearly 
exhibits nodal behavior (b) and the penetration depth stays 
linear down to low T in agreement with experiments (a, inset). 


spin-fluctuation mediated pairing. 

Spin-fluctuation pairing. What is the dominant pair¬ 
ing instability for the low-T orbitally ordered state? To 
answer this question, we consider the scattering vertex in 
the singlet channel projected onto the band space r(k, k') 
(see SI) and solve the linearized gap equation 



dS'T{k,k') 


ga(k') 

kFi(k')| 


Aaga(k), (5) 


where UFj(k') is the Fermi velocity of band j and the 
integration is performed over FSj to obtain the gap sym¬ 
metry functions pa (k) and the set of eigenvalues Aq . The 
largest eigenvalue corresponds to the leading instabil¬ 
ity and the corresponding eigenfunction determines the 
structure of the superconducting gap A(k) ^ ^(k) close 
to Tc- In order to solve Eq. <§ , the FS is discretized 
using a Delaunay triangulation pT] such that it reduces 
to solving a matrix eigenvalue problem. In the absence 
of orbital order, the leading instability is d-wave with 
nodes on the hole pockets, and no accidental nodes on 
the electron pockets, whereas in the absence of any band 
renormalization the leading instability is a nodeless sign 
changing s± state. In Fig.|^a) we show the result for the 
gap structure in the low-T phase with orbital order. The 
character of the gap structure cannot be classified in s or 
d-wave symmetry because the underlying band structure 
is only C 2 symmetric. [H] As seen Fig. Ba) the orbital 
order has strong effects on the position of the nodes, i.e. 
it removes the nodes from the hole pockets and induces 
nodal lines on the A-centered {Y for the other twin) elec¬ 
tron pocket. The associated DOS (maximum gap set to 
« 2.2 meV) and the linear-T behavior of the low-T pen¬ 
etration depth A shown in Fig. are remarkably similar 
to recent experimental findings. mm 

In summary, we have presented a model for the elec¬ 
tronic structure of FeSe that includes orbital ordering, 
which is consistent with recent ARPES and QO exper¬ 
iments on high quality FeSe samples. This band, along 
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with the standard local interaction potentials and ex¬ 
changes, explains both the T dependence of the NMR 
Knight shifts and spin-relaxation rate, and leads to a 
pairing state with nodes and a T dependence of the Lon¬ 
don penetration depth in agreement with a series of re¬ 
cent experiments on FeSe. 
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SUPPLEMENTARY INFORMATION = 0, 


Tight-Binding Bandstructure. Here we provide addi¬ 
tional information regarding the tight-binding parame¬ 
ters used for constructing the ten orbital model for FeSe. 
The tight-binding Hamiltonian for FeSe can be written 
as 

= £3 + 2 ^ 33(003 fci - 1 - cos fc 2 ) 

- 1 - 2^33 COs( 2 fc 2 ,) -1- 2<33 COs(2fcy) -1- 4^33 C0s(2fca;) cos(2fcy) 
-|-[ 2 t 33 ^ -1- 4 ^ 33 ^ COs(2fca;) -|- COs(2fcy)] COSfcz, 

Htb = V(k)4^(k)c^<T(k). (SI) 

1^,(7 

= 2 < 34 (cOsfci — COSfc 2 ), 

The orbitals /i, v follow the same notation convention 
used in Eschrig et al. In particular, the ten orbitals 
are [(xy)+, {ixz)+, iiyz)+, {z'^)+, [xy)-, (x^- 

y^)~, {—ixz)~, {—iyz)~, {z^)~], where the -I- and - refer 
to the two iron atoms in the unit cell. To reproduce 
the experimental results we renormalize the Hamilto¬ 
nian given in ReffHI by a factor of 6 and include addi¬ 
tional band shifts that lead to similar Fermi energies and 
Fermi pockets with shape and sizes similar to results of 
ARPES and quantum oscillation experiments. In general 
the Eermi surfaces seen in ARPES and QO experiments 
are best reproduced by including an additional spin or¬ 
bit coupling of Aoo ^ 20meV as shown in Eig.l of the 
manuscript. The corresponding case without the SO cou¬ 
pling is shown in Eig|V| below. One important difference 
between the two cases is the hole pockets at the zone cen¬ 
ter shown in Eig|V] become an outer hole cylinder when 
the effect of SO coupling is included in the tight binding 
model. Below we list the band dispersion given in RefflTI 
(with minor printing errors corrected). 

= 2itll{sinki -j- sinfc 2 ), 

= £3 + 2t33(cOS fci -1- COS fc 2 ) -1- 2t33 COs(2fca;) 
- 1 - 2^33 COs(2fcy) -1- 4^33 COs(2fca;) cos(2fcy), 

-|-[ 2 t 33 ^ -1- COs(2fca;) -1- 4 ^ 33 ^ COs(2fcy)] cos kz, 

= 2itll{sinki — sinfc 2 ), 

= £5 

H^q~ = 2 t}g(cos kx + cos ky) 

-l- 2 tig [(cos fci - 1 - cos fc 2 ) (cos kx + cos ky ) 

— sin fci (sin kx + sin fcj, ) -I- sin k 2 (sin kx — sin ky)] 
-l-4t[g^ (cos kx + cos ky) cos fcz 
-|- 2 t)g^{[cos(fci - 1 - ky) cos(fci - 1 - kx)] exp(ifcz) 
-l-[cos(fc 2 - 1 - ky) cos(fc 2 - kx)] exp(-ifcz)} 

= ei + (cos fci - 1 - cos ^ 2 ) 

-\-2tfl[cos{2ki) cos( 2 A: 2 )] 

-\-[2tii^ -|-4t}}^(cos ki -\-c 0 sk 2 ) 

-\-4:t1i^{cos2kx + cos2ky)] cos kz, 

H+- = 0 

= 2 it}g sin kx — 4(t}8^ sin kx + 4 °^ sin ky) sin kz 
-|- 2 jt[sin(fci - 1 - ky) exp(ifcz) — sin(fc 2 -f ky) exp(—zfcz)], 

= 0 , 

= 2itlg sin ky — 4(t[g^ sin kx + tig^ sin ky ) sin kz 

= 2 zt) 3 (sin ki — sin ^ 2 ) — sin( 2 A:j,) sin kz, 

-|-2it[g^[sin(fci -1- kx) exp(ifcz) -I- sin(fc 2 - kx) exp(-ifcz)], 

= 2 ztJ 3 (sinfci -|-sinA: 2 ) — 4^^°^ sin( 2 fca;) sin/cz, 

Htio = 0 

= 2^35 (cos fci — cosA: 2 ), 

= 2t^{coskx + cosky). 

= €2 + 2 t 22 (cOsfci -|- COSfc 2 ), 

= — 2^29 sin fcy, 

= 21^23 (sin fci + sinfc 2 ), 

= 2itlgsmkx, 

= 2it^{— sin fci -|- sin k 2 ), 

^ 2,10 = 2 t 2 %(cos kx - cos ky), 
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FIG. V. (Color online) Fermi surface without SO coupling for 
T > Ts (a) and T < Ts (b). As in Fig. 1 of the main text, 
the Fermi surfaces are derived from the ten orbital model, 
plotted in the crystallographic BZ and visualized with the 
summed-color method where the absolute value of the overlap 
is mapped to the RGB value of the color on the surface as 
indicated by the color bar. 

= 2 l 3 g cos kx + 2 i 4 g cos ky 

+213! [(cos ki + cos/cg) coskx — (sinfci — sinA;2) sinfca,] 
+2i4g[(cosfci + cos fcg) cosky — (sinfci + sinfc2) sinfcy] 
+4(f3g^ cos kx + cos ky ) cos k^ 
+ 2 i 3 |^[cos(fci + kx) exp{ikz) + cos(/c 2 — kx) exp(—f/cz)] 
+2<4g^[cos(fci + ky) exp(ffcz) + cos(fc2 + ky) exp(-jA;z)], 

^39~ ={cos kx + cos ky) sin kz, 

H3J0 = 2f4%sinfcy, 

= 2tlg cos{kx) + 2^8 cos{ky) 

+2149 [(cos ^1 + cos ftg) coskx — (sinfci — sinA;2) sinfca,] 
+2i3g[(cosfci + cos fcg) cosky — (sinfci + sinfc2) sinfc^] 
+ 4 (i 4 g^ cos kx + ^38^ COS ky ) cos fcz 
+ 2 i 4 g^[cos(fci + kx) exp(zfcz) + cos(fc 2 - kx) exp(-jfcz)] 
+2<3i^[cos(fci + ky) exp(jfcz) + cos(fc2 + ky) exp(-zfcz)], 

774 ig = 2 ft 4 ig sin fc^,, 

^5,10 — 

Below we provide the tight-binding hoppings that fit the 
DFT band structure as given in Eschrig et a7[T4] The 
hopping shifts performed in this study have been written 
as an added term to the original hopping integrals that 
undergo a renormalization by a factor oi z = 6. For the 
2D dispersion the hopping integrals are 

t\l = 0.086/z,<[g = -0.063/0 - 0.0211, 
tl° = -0.028/z + 0.0028, = 0.011/z, 

t\\ = -O.ObOi/z, til = 0.3051/0 -f 0.0761, 



FIG. VI. (Color online) Comparison of the band structure 
from the tight-binding model from ten orbitals (black line) 
and five orbitals (red line) (for T/ssand = -f ) show¬ 
ing the nearly perfect agreement in the high temperature 
tetragonal phase (a) and the low temperature orbitally or¬ 
dered state (b). 

til = - 0 . 109 / 0 , til = - 0 . 412/0 -f 0.022, 
til = -O.O66/0 - 0.003, til = -0.3641/0 - 0.0341, 

123 = 0.0891/0, tl% = 0.338/0 - 0.018, 
til = 0.232/0, til = O.O8O/0 -h 0.0025, 

^20 = 0 . 000lz,tll = O.OI6/0, 

= -0.045/0, tfg = 0.311/0 - 0.002, 

1^1 = 0 . 021 / z,til = -0.019/0, 
til = 0.099/0,4% = O.I8OI/0, 
til = O.146Z/0, 

Cl = 0.014/0, 
eg = -0.539/0-h 0.029, 

£3 = 0.020/0, 

£4 = 0.020/0, 

eg = -0.581/0-h 0.032 

For the fcz dispersion the additional hopping terms are, 

= 0 -h 0.0027, 

44 =0,44 = -0.017/0, 

411 = 0,44 = 0.0091/0, 

44 = 0.017/0 - 0.0027,44 = 0.0201/0, 

44 = 0.00301/0 -h 0.00451,44 = -0.00311/0, 



















tOOl = 0.011/z, ^8^ = 0 . 006 /z + 0.002, 

^201 = - 0 . 008 /z,i^^i = - 0 . 003 /z, 
tlf = 0 . 020 /z, 4 gi = 0 . 015 /z, 

= 0 . 025 /z - 0 . 0020 , 

421 = 0 . 006 /z. 

The five orbital tight-binding model can be generated 
from the ten orbital parameters given above by taking 
the symmetric and antisymmetric combinations as given 
in Refm The five orbital model is given by 

= ±H+-. (S 2 ) 

The tight-binding parameters of the five orbital 3 D model 
for the real orbitals {xy,x^ — ^xz,yz,z^) are given by 

t\\ = 0 . 086 /z,t}^ = - 0 . 063 /z - 0 . 0211 , 

^20 = _o. 028 /z -f 0 . 0028 , tH = 0 . 017 /z, 
t\l = - 0 . 056 i/z, = 0 . 305 i/z -f 0 . 076 i, 

til = - 0 . 109 /^, 4 ? = - 0 . 412 /z -f 0 . 022 , 
til = - 0 . 066 /z - 0 . 003 , = - 0 . 364 i/z - 0 . 034 i, 

til = 0-089i/z,tl% = 0 . 338 /z- 0 . 018 , 
til = 0.282/ z, til = 0 . 080/0 -f 0 . 0025 , 

^20 = 0 . 009 / 0,^21 = 0 . 016 / 0 , 
til = -0.045/0,= 0.311/0 - 0.002, 
til = 0.027/z,tll = -0.019/0, 

41 = 0.099/0,4% = -O.18Oi/0, 

^5 = -O.146i/0, 
ei = 0 . 014 / 0 , 

£2 = -0.539/0-f 0 . 029 , 

£3 = 0.020/0, 

£4 = 0.020/0, 

£5 = - 0 . 581 / 0 -f 0 . 032 . 


For the dispersion the additional hopping terms are 


= 0.0-f 0.0165, 

44=0.0,44 = 0.0, 

44 = 0.0,44 = 0.0, 

44 = 0 . 017/0 - 0.0165,44 = 0 . 020 / 0 , 
44 = 0.0030/0 -f 0.027,44 = 0.0, 
44 = 0 . 011 / 0,44 = 0 . 006 / 0 -f 0.0152, 
44 = -0.008/0,44 = 0.0, 

44 = 0.0,44 = 0 . 015 / 0 , 

44 = 0 . 025/0 - 0 . 0122 , 

44 = 0.0. 


The above five orbital fit agrees well with the ten or¬ 
bital model as can be seen from Fig|VI| In Fig|VII| we 
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FIG. VII. (Color online) Extremal orbit area of Fermi pocket 
at low temperatures in the absence of SO coupling. The plots 
are for various 8 values and shows the small difference in the 
dispersion of the Fermi pockets between cj) = 0° (a) and 
4> = 90° (b). Here 6 is the angle with c-axis and <(> = 0° is 
aligned along the Fe-Fe direction. 


show the QO results for magnetic fields along two or¬ 
thogonal in-plane directions. As can be seen in Fig |VII| 
the fcj dispersion [9 dependence) of the orbit areas are 
slightly different for </ = 0° and </> = 90°. These two 
(j) values would correspond to the areas observed from 
the two structural domains for a given magnetic field di¬ 
rection in quantum oscillation experiment and therefore 
have slightly different QO frequencies when 9 jtz Q, 

Interaction Hamiltonian. The local interactions are 
included via the Hubbard-Hund Hamiltonian 

H = Hq + U'^ nu^nai + U' nun^t 

4£cr4£'cr'Gfcr'Ci^V (S3) 

i,i' <l <y,a' 

where the interaction parameters U, 11/ J, J' are given 
in the notation of Kuroki et al. m- Here, I is an orbital 
index with £ S (1,..., 5) corresponding to the Fe-orbitals. 
The pairing interaction in band representation is now cal¬ 
culated from the interaction in the orbital representation 
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FIG. VIII. (Color online) RPA susceptibility in the tetragonal 
phase at r = 110 K calculated for Uz = 1.8 eV and Jz = 
Q.lUz and plotted in the Brillonin zone at fcz = 0 (a) and 
along a high symmetry path in comparison with the non¬ 
interacting susceptibility (b). 


via 


r(k,k')=Re ^ a5(’*(k)al:’*(-k) (S4) 

x[r,,,,,3,,(k,k')] af^.(k')al5(-k'), 

where k and k' are quasiparticle momenta restricted to 
the pockets k G FSi and k' G FSj, where i and j corre¬ 
spond to the band index of the Fermi surface sheets. The 
vertex function in orbital space describes the 

particle-particle scattering of electrons in orbitals £ 2,^3 
into £ 1 ,^ 4 . 

In RPA it is given by 








■ (S5) 




Note further, that in the spin-singlet channel, we 
symmetrize the pairing vertex by using the expres¬ 
sion l/2[F(k,k') - 1 - r(k, —k')] in the linearized gap 
equation[5S]. The RPA susceptibility in q-space is shown 


FIG. IX. (Color online) RPA susceptibility in the orbital or¬ 
dered phase at T = 40K (b) calculated for Uz = 1.8 eV and 
Jz = O.lUz and plotted in the Brillouin zone at fcz = 0 (a) 
and along a high symmetry path in comparison with the non¬ 
interacting susceptibility (b). 




FIG. X. (Color online) Superconducting gap function plotted 
over the Fermi surface of the 5 orbital model for Uz = 1.92 eV 
and Jz = 0.125Uz (a) and corresponding density of states (b). 


in Figs IVIII] and The w = 0 susceptibility is domi¬ 
nated by the (tt, tt) nesting between the electron Fermi 
pockets. Note, however, that due to the small band en¬ 
ergies for this system, the q-dependence at relatively low 
finite uj can be quite different, as we will explore else¬ 
where. The solution of the linearized gap equation for 
Uz = 1.92 eV and Jz = O.125C/0 is plotted in Fig|^ 
It can be seen that in comparison to the gap struc¬ 
ture obtained for Uz = 1.8eV and Jz = Q.lUz in the 
manuscript, the gap structure for a larger U leads to a 
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FIG. XL (Color online) Spin lattice relaxation rate calculated 
in the absence of structure factor filtering i.e Ah/(q) = 1 for 
the same parameters as used in Fig.3 of the manuscript. 


weakening of the nodes on the electron pocket. The den¬ 
sities of states p{ijj) are then calculated using the gap 
A(k) = Ao 5 (k) interpolated on a fine fc-mesh of approx¬ 
imately 1000 x 1000 x 200 points and an imaginary smear¬ 
ing of 77 = 0.02 meV. 

Finally, the effect of the density of states on the pene¬ 
tration depth is A is estimated by the formula 

A oc y dojp{ijj) y ^ , (S 6 ) 


where /(w) is the Fermi function. 

Spin Lattice Relaxation Rate. The NMR spin lattice 
relaxation rate is given as 


1 

Tif 


1- / M2lin[XfiPA(q,w)] 

a ^ Ahf{q) - - 

(.j—>0 ^' ui 


(S7) 


where XflPA(q, w) is the RPA enhanced susceptibility 
and Ahf{q) = A^j cos{qoc) cos{qy) is the hyperfine form 
factor which is a 3 x 3 matrix with non zero diagonal 
elements in the paramagnetic state. The q dependence 
in the form factor leads to filtering of the fluctuations at 
the Brillouin zone edges. Comparing the relaxation rate 
given in Fig.3(b) in the manuscript with the correspond¬ 
ing case without the form factor (see Fig 0 we find that 
in the absence of filtering for Ahf{c[) = 1 , the upturn in 
l/TiT takes place at a higher temperature. 

For the multi-orbital model, the expression for 1/TiT 
can be simplified by expanding the bare susceptibility 
in powers of w. After some algebra, this leads to the 
following expression. 


1 

7\T 


mat.el. q 


^(q) 

i - f7i(q) 


-UBia) ' 


1 - UA{q) 


1 - UA{q) 


(S 8 ) 


where the matrices A, and B are given by 
^(q) = X ^M.^(k,q)R-e(G^.j.(k,q,0)) 

k./j.i/' 



T(K) 



FIG. XII. (Golor online) Plot showing the effects of orbital or¬ 
dering on spin-lattice relaxation rate and NMR Knight shift. 
The hyperfine form factor Ahf = a -I- P{T) are same as in 
manuscript Fig. 3 except the temperature dependent factor 
P{T) = 0 in the absence of orbital order, (a) NMR Knight 
shift versus T comparing cases with and without orbital or¬ 
der. (b) Spin-lattice relaxation rate versus T for cases with 
and without orbital order. 


B(q) 


x[/(F;.(k-t-q))-/(F;^(k))], 

x[/(F;.(k-|-q))-/(F;^(k))]. 


duj 


(S9) 

u;—0 

(SIO) 


Here q) = a^(k)aP*(k)a« (k + q)a^*(k -|- q) is in 

general complex with a® (k) representing the eigenvector 
component corresponding to band p and orbital s, 
is the eigenvalue for band p, and f{Ey(k)) is the Fermi 
function at energy iil^(k). The Greens function reads 


G^,i.(k,q,a;) = 


E^(k + q) 


1 

Ef^(k) — uj — iS 


(Sll) 


The above expressions for Knight shift and spin lattice 
relaxation rate have been evaluated for our band struc¬ 
ture in Fig.3 of the manuscript. The effect of orbital 
ordering can be made clearer upon comparing the cases 
with and without orbital ordering. As shown in Fig |XII| 
(a), we find that in the absence of orbital ordering the 
NMR Knight shift saturates as we go to low temper¬ 
atures. The spin lattice relaxation rate in Fig XII (b) 
shows that although band structure effects do cause a 
weak enhancement of spin fluctuations even in the ab¬ 
sence of orbital ordering, the effect is enhanced due to 
the presence of orbital ordering. 

s+d wave orbital order. 

The effect of an additional d-wave orbital order can be 
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FIG. XIII. (Color online) (a) Band structure of the 10 orbital 
model at low T with d-wa,Ye orbital order and SO coupling, 
yielding the QO frequencies as a function of magnetic field an¬ 
gle 9 shown in (b), where the error bars indicate the numerical 
uncertainty in the determination of the extremal orbits. The 
orbital character in the band plot is indicated by the colors 
red dxz, green dyz, blue d^y, yellow d^ 2 _yi, purple d^^ 2 _^ 2 . 


maximum value of 25meV at T = 10 K). This leads to 
the hole band near T being of primarily dxz idyz) charac¬ 
ter, and the inner electron band at M point of mainly dyz 
(dxz) character as can be seen in the orbitally resolved 
band shown in Fig |XIIl| This opposite shift of the bands 
can be seen more clearly from Fig XIV (inset) which dis¬ 
plays the difference in band energy {AE = Exz{T = 
170K)) — Exz{T = 10K)) between the dxz band in or¬ 
bital ordered state at T = 10 K) and no orbital ordering 
at T = 170 K). It shows clearly that the band shift due 
to orbital ordering moves in opposite directions at the 
F and X points, and therefore undergoes a sign change 
along the F—X direction (in the unfolded Brillouin zone). 
Unlike the Fermi surface of the case with pure s-wave or¬ 
bital order (Fig. |^b)), the opposite shift of the band 
(at F and X points) due to the s + d-wa.ve orbital or¬ 
der leads to different orbital content between the hole 
and the inner electron Fermi surface pockets, as seen in 
the Fermi surface plot in Fig |XIV| (b). The extremal or¬ 
bits measurable by QO have also been calculated in the 
presence of d-wave orbital order. As seen in Fig |XIIl[ 
the Fermi surface areas agree well with experiments sim¬ 
ilar to previous results for purely s-wave orbital order. 
Additionally the Sommerfeld co-efficient calculated from 
the effective masses derived from the quantum oscilation 
calculation for d-wave orbital order gives a value of 4.1 
mJ/mol-K^. Calculations of the NMR Knight shift and 
spin-lattice relaxation rate 1/TiT in the presence of the 
d-wave orbitally ordered state are shown in Fig |XV| Evi¬ 
dently, general agreement with bulk experiments can also 
be obtained for this case. 


included with the term 

Hoo = y^[As(r)(nxzer(k) - Hyza-Oi-)) + 

kcr 

Ad{T){cos{kx) - cos{ky)){nxza{^) + ny2<^(k))],(S12) 

where both As(T) and Ad(T) exhibit a typical mean-field 
T dependence and have been chosen such that Ay/Ad > 
0. The band structure shifts presented above have to be 
slightly adjusted to agree with this orbital order. The 
additional band shifts are 

Atlf = At°“ = 0.00016, 

Atll = Atf^ = -0.0013, 

At}° = -Atll = -0.0006. 

For As{T) > 0 and Ad{T) > 0, the s- and d-wave 
orbital components move the dxz/dyz bands in opposite 
directions at the F and M points, respectively. Both 
As{T) and Ad{T) follow a mean field behavior with a 



FIG. XIV. (Golor online) Fermi surface for T > Ts (a) and 
T < Ts (b) in the presence of a d-wave orbital order. As in 
Fig. 1 of the main text, the Fermi surfaces are derived from 
the 10 orbital model, plotted in the crystallographic BZ and 
visualized with the summed-color method where the absolute 
value of the overlap is mapped to the RGB value of the color 
on the surface as indicated by the color bar. The inset shows 
the shift in energy of the d^z band between T=170K and 
r=10K due to the effect of a s-|-d-wave orbital order. The d^z 
band (and similarly the dyz band) shifts in opposite directions 
at the r and A-points. 
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FIG. XV. (Color online) Plot showing the effects of d-wave or¬ 
bital ordering on spin-lattice relaxation rate and NMR Knight 
shift, (a) NMR Knight shift versus T. (b) Spin-lattice relax¬ 
ation rate versus T. 
















